/*
 * topological-hex, a program to compute combinatorial hexahedral meshes.
 *
 * Copyright (C) <2018> <Université catholique de Louvain (UCL), Belgique>
 *
 * List of the contributors to the development, description and complete
 * License: see LICENSE file.
 *
 * This program (topological-hex) is free software:
 * you can redistribute it and/or modify it under the terms
 * of the GNU General Public License as published by the Free
 * Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program (see COPYING file).  If not,
 * see <http://www.gnu.org/licenses/>.
 */


#include <check.h>
#include <stdlib.h>
#include <stdbool.h>

#include "mesh.h"
#include "input.h"

static bool find_quad(const mesh *mesh, const vertex_index *quad) {
  for (quad_index i = 0; i < mesh->num_boundary_quads; i++) {
    quad_vertex_index j0;
    for (j0 = 0; j0 < 4 && quad[j0] != mesh->quads[4*i]; j0++);
    if (j0 == 4)
      continue;

    quad_vertex_index j1 = (j0 + 1) % 4;
    quad_vertex_index j2 = (j0 + 2) % 4;
    quad_vertex_index j3 = (j0 + 3) % 4;

    if (quad[j0] == mesh->quads[4*i + 0] &&
        quad[j1] == mesh->quads[4*i + 1] &&
        quad[j2] == mesh->quads[4*i + 2] &&
        quad[j3] == mesh->quads[4*i + 3])
      return true;

    if (quad[j0] == mesh->quads[4*i + 0] &&
        quad[j3] == mesh->quads[4*i + 1] &&
        quad[j2] == mesh->quads[4*i + 2] &&
        quad[j1] == mesh->quads[4*i + 3])
      return true;
  }

  return false;
}

static void check_quads(const mesh *mesh, uint32_t n,
                        const vertex_index (*quad)[4]) {
  for (uint32_t i = 0; i < n; i++) {
    ck_assert_msg(find_quad(mesh, quad[i]),
                  "Missing quad: (%u, %u, %u, %u)",
                  quad[i][0], quad[i][1], quad[i][2], quad[i][3]);
  }
  ck_assert_int_eq(mesh->num_boundary_quads, n);
  ck_assert_int_eq(mesh->num_quads, n);
}

START_TEST(test_load_error) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/does_not_exist.mesh") == IO_ERROR);

  mesh_release(&mesh);
} END_TEST

START_TEST(test_load_tetrahedron) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/tetrahedron.mesh") == SUCCESS);

  ck_assert_int_eq(mesh.num_vertices, 14);
  const vertex_index quads[][4] = {
    {0, 7, 4, 9},
    {1, 8, 4, 9},
    {2, 7, 4, 8},
    {1, 8, 6, 12},
    {3, 12, 6, 13},
    {2, 8, 6, 13},
    {0, 7, 5, 10},
    {2, 7, 5, 13},
    {3, 10, 5, 13},
    {0, 9, 11, 10},
    {1, 9, 11, 12},
    {3, 10, 11, 12},
  };
  check_quads(&mesh, sizeof(quads)/sizeof(*quads), quads);

  mesh_release(&mesh);
} END_TEST

START_TEST(test_load_hexahedron) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/hex.mesh") == SUCCESS);

  ck_assert_int_eq(mesh.num_vertices, 26);
  const vertex_index quads[][4] = {
    {0, 8, 20, 16},
    {1, 8, 20, 18},
    {3, 11, 20, 18},
    {2, 11, 20, 16},
    {0, 9, 21, 16},
    {4, 9, 21, 17},
    {2, 12, 21, 16},
    {6, 12, 21, 17},
    {1, 10, 22, 18},
    {5, 10, 22, 19},
    {7, 13, 22, 19},
    {3, 13, 22, 18},
    {0, 8, 23, 9},
    {1, 8, 23, 10},
    {5, 10, 23, 14},
    {4, 9, 23, 14},
    {2, 11, 24, 12},
    {3, 11, 24, 13},
    {7, 13, 24, 15},
    {6, 12, 24, 15},
    {6, 15, 25, 17},
    {7, 15, 25, 19},
    {5, 14, 25, 19},
    {4, 14, 25, 17},
  };
  check_quads(&mesh, sizeof(quads)/sizeof(*quads), quads);

  mesh_release(&mesh);
} END_TEST

START_TEST(test_load_hex_mesh) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/schneiders_88.mesh") == SUCCESS);

  ck_assert_int_eq(mesh.num_vertices, 105);
  ck_assert_int_eq(mesh.num_hexes, 88);

  mesh_release(&mesh);
} END_TEST

START_TEST(test_reserve) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/hex.mesh") == SUCCESS);

  ck_assert_int_eq(mesh.num_boundary_quads, 24);
  vertex_index quad_copy[24*4];
  memcpy(quad_copy, mesh.quads, sizeof(quad_copy));

  ck_assert(mesh_reserve(&mesh, 79) == SUCCESS);

  ck_assert(memcmp(quad_copy, mesh.quads, sizeof(quad_copy)) == 0);
  ck_assert_int_eq(mesh.num_boundary_quads, 24);
  ck_assert_int_eq(mesh.num_quads, 249);
  ck_assert_int_eq(mesh.num_hexes, 79);

  for (uint32_t i = 4*mesh.num_boundary_quads; i < 4*mesh.num_quads; i++)
    mesh.quads[i] = 5;

  for (uint32_t i = 0; i < 8*mesh.num_hexes; i++)
    mesh.hexes[i] = 3;

  for (uint32_t i = 0; i < 6*mesh.num_hexes; i++) {
    mesh.hex_adjacency[i] = 4;
    mesh.hex_facets[i] = 7;
  }

  mesh_release(&mesh);
} END_TEST

START_TEST(test_reserve_smaller_mesh) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/schneiders_88.mesh") == SUCCESS);
  ck_assert(mesh_reserve(&mesh, 64) == SUCCESS);
  mesh_release(&mesh);
} END_TEST

START_TEST(test_bad_inputs) {
  mesh mesh;
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/mismatch-hex-count.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/mismatch-quad-count.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/mismatch-vertex-count.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/vertex-negative.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/vertex-too-large.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/no-hex-count.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/no-quad-count.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/no-vertex-count.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/not-enough-hexes.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/not-enough-quads.mesh") != SUCCESS);
  ck_assert(mesh_load(&mesh, INPUT_DIR "/malformed_data/not-enough-vertices.mesh") != SUCCESS);
} END_TEST

Suite *mesh_suite(void) {
  Suite *s = suite_create("Mesh");

  TCase *tcase = tcase_create("mesh_load");
  tcase_add_test(tcase, test_load_error);
  tcase_add_test(tcase, test_load_tetrahedron);
  tcase_add_test(tcase, test_load_hexahedron);
  tcase_add_test(tcase, test_load_hex_mesh);
  suite_add_tcase(s, tcase);

  tcase = tcase_create("mesh_reserve");
  tcase_add_test(tcase, test_reserve);
  tcase_add_test(tcase, test_reserve_smaller_mesh);
  suite_add_tcase(s, tcase);

  tcase = tcase_create("bad inputs");
  tcase_add_test(tcase, test_bad_inputs);
  suite_add_tcase(s, tcase);

  return s;
}
